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Key words Abstract Fungi are rich in complexes of cryptic species that need a combination of different approaches to be 
Beauveria delimited, including genomic information. Beauveria (Cordycipitaceae, Hypocreales) is a well-known genus of 
chemotaxonomy entomopathogenic fungi, used as a biocontrol agent. In this study we present a polyphasic taxonomy regarding 
population genomics two widely distributed complexes of Beauveria: B. asiatica and B. bassiana s.lat. Some of the genetic groups as 
taxonomy previously detected within both taxa were either confirmed or fused using population genomics. High levels of 


divergence were found between two clades in B. asiatica and among three clades in B. bassiana, supporting their 
subdivision as distinct species. Morphological examination focusing on the width and the length of phialides and 
conidia showed no difference among the clades within B. bassiana while conidial length was significantly different 
among clades within B. asiatica. The secondary metabolite profiles obtained by liquid chromatography-mass spectro- 
metry (LC-MS) allowed a distinction between B. asiatica and B. bassiana, but not between the clades therein. Based 
on these genomic, morphological, chemical data, we proposed a clade of B. asiatica as a new species, named 
B. thailandica, and two clades of B. bassiana to respectively represent B. namnaoensis and B. neobassiana spp. 
nov. Such closely related but divergent species with different host ranges have potential to elucidate the evolution 
of host specificity, with potential biocontrol application. 
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INTRODUCTION 


Fungi are rich in complexes of cryptic species, i.e., morphologi- 
cally similar but genetically isolated. In pathogenic fungi, cryptic 
sibling species are often specialized on different hosts (Le Gac 
et al. 2007, Giraud et al. 2010, Kobmoo et al. 2012), so that it 
is essential to elucidate their genetic subdivision and species 
limits. Taxonomy is therefore an essential tool for understanding 
the evolution of host specificity in pathogens, and to evaluate 
the risk of attacking untargeted species with biocontrol agents. 


Beauveria bassiana was first discovered in the 18th century as 
representative of the genus Botrytis (Basalmo-Crivelli 1835), 
before being renamed as the first species of the genus Beau- 
veria (Vuillemin 1912). It is distributed worldwide (Rehner & 
Buckley 2005, Imoulan et al. 2017). Being recognised as a rich 
source of efficient mycoinsecticides, in particular B. bassiana 
and B. brongniartii (Zimmermann 2007), the genus Beauveria 
has received a lot of attention as potential biological control 
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agents (García-Estrada et al. 2016). The discovery of new 
species and their precise delimitation contribute to enlarge 
the potential for finding suitable biocontrol agents and for their 
safe use. Beauveria is characterised by the formation of short, 
sympodial, globose to flask-shaped phialides with holoblastic 
conidia (Khonsanit et al. 2020). The size and the shape of co- 
nidia are variable within the genus and historically constituted 
discriminant characters for delimitating and identifying species 
(Rehner & Buckley 2005, Imoulan et al. 2017, Abdessamad 
2019). However, these morphological characters can exhibit 
overlapping values between species, even between relatively 
distantly related ones (Khonsanit et al. 2020). Furthermore, 
many new Beauveria species have been proposed in recent 
years, mainly based on the monophyly of clades in molecular 
phylogenies, but sometimes with only a few samples, without 
clear morphological distinction (Chen et al. 2018). Beauveria is 
thus a genus with potential cryptic species waiting to be eluci- 
dated. However, there has never been a thorough morphometric 
analysis with solid statistical tests. 


Although species of Beauveria are known to produce an array 
of secondary metabolites (Xu et al. 2008, 2009, Rohlfs & 
Churchill 2011, Udompaisarn et al. 2020), chemical compounds 
produced have been rarely used for taxonomic purpose in this 
genus; the potential of using chemical diversity as markers for 
classification and species identification was explored several 
years ago for Beauveria (see Mugnai et al. 1989, Bridge et 
al. 1990) without any update for recently discovered species. 
Chemical compounds allowed the distinction between popula- 
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tions within B. bassiana (Bridge et al. 1990) but this finding was 
not compared with genetic data. 


The abundant genomic resources presently available for Beau- 
veria (Valero-Jiménez et al. 2016, Toopaang et al. 2017, Mei 
et al. 2020), combined with decreasing cost for whole-genome 
sequencing (WGS), make it possible to access novel genetic 
information on a large number of samples with reasonable 
cost and time. Next generation sequencing (NGS) has in fact 
recently been used to reveal fungal cryptic species (Sepulveda 
et al. 2017, Kobmoo et al. 2019, Matute & Sepúlveda 2019). 


In the present study, we aimed at revising the taxonomy of two 
important Beauveria taxa, B. bassiana and B. asiatica. In the 
past, populations of B. bassiana were studied using a limited 
numbers of markers (Wang et al. 2003, 2005, Mitina et al. 
2011, Meyling et al. 2012), some of which are now considered 
unreliable, such as Amplified Fragment Length Polymorphism 
(AFLP) (Aquino De Muro et al. 2005, Trissi et al. 2013). A re- 
cent population genomics study of these species showed clear 
intraspecific genetic groups, but their taxonomic status was not 
evaluated with appropriate phylogenetic frameworks or detailed 
examination of their morphology and chemistry (Mei et al. 2020). 


In our previous taxonomic revision of Beauveria, we found 
poorly supported monophyletic clades based on multigene 
phylogenies within B. asiatica (four clades named A, B, C and 
D) and B. bassiana (three clades named A, B and C) (Khonsanit 
et al. 2020). The poor node supports did not allow the recog- 
nition of these clades as distinct species, and morphological 
traits appeared overlapping. The poor node support may have 
been due to the limited number of markers. In the current study, 
our objectives were to evaluate the taxonomic status of the 
genetic groups within B. asiatica and B. bassiana by a multi- 
disciplinary approach including whole-genome sequencing 
(WGS), as well as statistical analyses of morphological charac- 
ters (the width and the length of conidia and phialides) and 
chemical profiles determined by liquid chromatography-mass 
spectrometry (LC-MS). A thorough taxonomic revision should 
indeed be based on a combined dataset including sufficient 
genetic information, as well as morphological and chemical 
traits on a high number of individuals. This constitutes the 
foundation of integrative (or holistic) taxonomy of fungi. 


MATERIALS AND METHODS 


Fungal culture and DNA extraction 


Seventy-eight strains of presumably 12 Beauveria species from 
the BIOTEC culture collection (BCC) and ARS oollection of 
entomopathogenic fungal cultures (ARSEF), including the ex- 
type strains of known species and strains previously identified as 
B. asiatica or B. bassiana as well as some unidentified strains fol- 
lowing Khonsanit et al. (2020), were used in this study (Table S1). 
The strains were cultured in 50 mL of potato dextrose broth and 
incubated at 25 °C for a week. Fungal mycelia were harvested 
by filtering with a sterilised nylon mesh and washed with ethyl- 
ene diamine tetraacetic acid (EDTA) and distilled water. DNA ex- 
traction was done using a cetrimonium bromide (CTAB)-based 
method following Kobmoo et al. (2019); the DNAs were purified 
using the high pure PCR template preparation kit (Roche). The 
quality and quantity of DNA were verified using Nanodrop™ One 
Microvolume UV-Vis Spectrophotometry (Thermo Fisher) and 
an electrophoresis on 0.8 % agarose gel at 100 V for 1/2 h. 


DNA library construction and sequencing 


For whole genome shotgun sequencing of 78 strains, approxi- 
mately 300 ng of each DNA sample was used for a library con- 
struction following the protocol in the MGIEazy FS library prep 
kit (MGI Tech, Shenzhen, China). The samples were pooled and 
sequenced in a single lane. Paired-end (150 bp) sequencing 


was performed on the MGISEQ-2000RS (MGI Tech, Shenzhen, 
China) according to the manufacturer's instructions. 


Data processing and detection of single nucleotide poly- 
morphisms (SNPs) 


Using the MegaBOLT v. 1.5.6.11 software package, raw reads 
were de-multiplexed according to their barcodes and the 
adapter/barcode sequences were removed. After removing 
the low-quality regions, clean reads were mapped to the 
B. bassiana reference genome ARSEF8028 (Valero-Jiménez 
etal. 2016) using the MegaBOLT v. 1.5.6.11 alignment software 
(Minimap2 v. 2.11-r797-v03) and the variants were called using 
GATK HaplotypeCaller v. 3.8. SNP markers with poor quality 
data were filtered out using the following criteria: 


— a minor allele frequency « 0.1; 
— depth coverage less than 10 x; 
— more than 20 % missing data. 


The ploidy was set as haploid as all the strains were collected 
from asexual mycelia, i.e., putatively haploid for this asco- 
mycete fungus. The raw reads data were deposited at NCBI 
Sequence Read Archive associated to the BioProject accession 
PRJNA744643. The SNPs data were deposited at Mendeley 
Data repository (Sonthirod et al. 2021). 


Population structures, genetic diversity and linkage dis- 
equilibrium 

The selected SNPs (729549 SNPs) were aligned and subjected 
to a maximum-likelihood based phylogenetic tree inference 
using RaxML v. 8 (Stamatakis 2014), by specifying GTRCAT 
model with Lewis biases ascertainment correction. The reli- 
ability of the resulting tree was evaluated with 1000 bootstraps. 
Bayesian clustering analyses were conducted with FastStruc- 
ture (Raj etal. 2014) and a principal component analysis (PCA) 
was achieved using the package adegenet in R (Jombart 2008). 
These analyses were intended to give an initial insight into 
interspecific divergence and intraspecific population structures, 
in order to assess the taxonomic status of the samples included 
as well and to reveal any genetic group that can potentially be 
identified as new species. We calculated the fixation index (F,,) 
and metrics of absolute divergence (D,,) between B. asiatica 
and B. bassiana, as well as between intraspecific genetic clus- 
ters as revealed by previous analyses. The nucleotide diversity 
(Pi) within species and genetic clusters were also calculated 
using the PopGenome package in R (Pfeifer et al. 2014). The 
r? (square of correlation coefficient representing statistical as- 
sociation between pairwise SNPs) among isolates belonging to 
B. asiatica, B. bassiana and their respective intraspecific genetic 
clusters were calculated between pairs of SNPs using VCFtools 
(Danecek et al. 2011). The decay of linkage disequilibrium was 
visualized using R. Neighbour-net phylogenetic network based 
on p-distance was constructed using the software SplitsTree 
v. 4.14.18 (Huson & Bryant 2006). 


Species tree 


Gene sequences were obtained using FasterAlternateRefe- 
renceMaker tool from GATK (McKenna et al. 2010) which 
altered the reference sequences with SNPs and simple indels; 
all complex substitutions were masked as Ns to reduce am- 
biguous alignment. The sequences were then aligned using 
the software MAFFT (Katoh & Standley 2013). A total of 1132 
genes were first selected based on: 

— their size (1 Kbp minimum and 10 Kbp maximum) to ensure 
an absence of recombination within genes; 

— their physical distance on the genome; only genes with at 
least 20 Kbp between the start and the end of consecutive 
genes were selected to avoid the linkage disequilibrium; 

— the number of SNPs per bp (» 0.02) to provide sufficient 
phylogenetic signal. 
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The alignments of the 1 132 genes were subjected to Bayesian 
analyses with MrBayes (Ronquist et al. 2012), each with 2 Markov 
chain Monte Carlo (MCMC) runs for 1 M generations with a 
20 96 burn-in phase, allowed to sample across the substitu- 
tion model space. Only the final gene trees with « 0.02 split 
frequencies were retained for 1034 genes. The 1034 genes 
sequences were also concatenated and subjected to MrBayes 
with the model GTRGAMMAI, as chosen by ModelTest-NG 
(Darriba et al. 2020), for 5 M of MCMC generations. The final 
tree was obtained after a 25 % burn-in phase. 


Bayesian Concordance Analyses (BCA) 


To evaluate whether the genetic clusters within species were 
consistently recovered by different markers throughout the ge- 
nome, Bayesian concordance analyses (BCA) were conducted 
for the three clusters of B. bassiana with ARSEF7032 (B. kipu- 
kae) as an outgroup, and for the two clusters of B. asiatica with 
ARSEF617 (B. brongniartii) as an outgroup. By focusing on 
each taxon of interest with its respective outgroup, the number 
of genes with variable SNPs naturally reduced. The BCA were 
thus done by selecting, among the 1 034 concatenated genes, 
only those which had at least two SNPs with no missing data 
and the number of SNPs per bp > 0.001, resulting in 887 genes 
for B. asiatica and 100 genes for B. bassiana. The single-gene 
trees obtained from Bayesian inferences were processed into 
BUCKy (Larget et al. 2010) which evaluated the concordance 
between the selected loci for the clades by giving genome-wide 
concordance factors with 95 96 interval which can be 0 (absence 
of concordance) to 1 (total concordance). 


Secondary metabolites profiling 


Strains of Beauveria, selected to represent the various genetic 
groups within B. asiatica and B. bassiana as well as ex-type 
strains of other species (Table S1), were grown in 200 mL of 
yeast with malt extract and glucose (YMG) medium (10 g of 
malt extract, 4 g of D-glucose, 4 g of yeast extract and 1000 
mL of distilled water, pH 6.3), incubated at 25 °C under shaking 
condition (140 rpm). Since growth rate varies among strains, the 
prolonged fermentation after glucose depletion for an individual 
strain was set to be half the time required for that particular 
strain to reach the glucose depletion (Table S2). The glucose 
content of each fermented broth was estimated using urine 
glucose test strips (DIRUI®, Jilin, China). The mycelia were then 
separated from the broth either by filtration or centrifugation. 


The extraction of fungal secondary metabolites was performed 
according to Phainuphong et al. (2017) with culture filtrates 
passing through successive extraction using acetone, ethyl 
acetate and methanol, resulting in a final methanol-based cell 
extract (CE). A single extract per strain was obtained. The ex- 
periment of liquid chromatography-mass spectrometry (LC-MS) 
was done using ultra-high performance liquid chromatography 
(UHPLC) - Orbitrap Fusion™ Tribrid™ mass spectrometer 
(Thermo Scientific, Massachusetts, USA), equipped with elec- 
trospray ionization (ESI) source. The separation of compounds 
was done on an Acquity UPLC? HSS T3 C18 column (1.8 um 
diam, 2.1 x 100 mm) maintained at 40 °C with flow rate of 
0.4 mL/min. A mobile phase system was solvent A (water with 
0.1 % formic acid) and solvent B (acetonitrile with 0.1 % for- 
mic acid). The optimized gradient for the best separation was 
0-12 min, 12-95 96 B; 12—14 min, 12—95 % B; 14—16 min, 
12-95 % B. The mass spectrometer was operated in ESI posi- 
tive which covered more than 80 96 of metabolites present in 
the tested pool sample. Each sample was injected in triplicate. 
The MS calibration was conducted using Pierce LTQ Velos ESI 
Positive lon Calibration (PSP3A 88323) and Pierce ESI Nega- 
tive lon Calibration (PSP3A 88324) according to manufacturer's 
protocol. Each extract was injected three times into the LC-MS 
as technical replicates. 


The software Compound Discoverer v. 3.1 was used for data 
pre-processing steps. The MS raw files were subjected to peak 
alignment, peak picking, adduct grouping and normalization 
with parameters adjusted to fit to the chromatographic data 
obtained from this experiment. The area under curve (AUC) of 
each metabolite was determined and normalized with that of 
pooled samples. The features (m/z at specific retention time 
(rt)) with relative standard deviations over 30 % and average 
group area in non-inoculated YMG media over 5 x 10? were 
filtered to ensure good-quality peaks and to cover peaks apart 
from culture media. Top 200 highest abundant features in the 
CE pool sample were selected for statistical analyses including 
a Euclidean distance-based neighbour-joining (NJ) tree and 
principal component analysis (PCA) using respectively the 
packages ‘ape’ (Paradis & Schliep 2019) and ‘FactoMineR’ (Lé 
et al. 2008) in R (R Core Team 2020). 


Conservation of secondary metabolites gene clusters 


We examined whether the distribution of secondary metabolites 
production among clades and species was reflected by the 
conservation pattern of secondary metabolites gene clusters 
(SMGCs). First, we inferred and annotated SMGCs present in 
the reference genome (ARSEF8028) using antiSMASH fungal 
v. 4.1.0 (Blin et al. 2017). The deduplicated and filtered mapped 
reads used in the phylogenomics above were proceeded to the 
pipeline of CNVnators (Abyzov et al. 2011) in order to detect 
regions of deletion in each isolate, based on normalised read 
depth. The detected blocks of deletion were intersected with 
the position of SMGCs and any SMGC with more than 50 % 
of its length overlapping a deletion is considered as absent or 
non-functional. The SMGCs presence/absence patterns were 
used to infer a cladogram representing the clustering of the 
isolates of Beauveria based on binary distances. 


Morphological examination and species description 


A thorough morphological investigation was conducted follow- 
ing the techniques described in Khonsanit et al. (2020). The 
samples and culture plates were photographed using a digital 
Nikon D5100 camera. Colony characteristics and microscopic 
measurements of phialides and conidia (length and width) were 
done after growth on potato dextrose agar medium (PDA, 20 g 
Difco potato dextrose agar, 1 L distilled water) and incubation 
under white light/dark cycles at room temperature; 30 phialides 
and conidia were measured for each strain. To statistically 
test the difference between the genetic groups within each 
of B. asiatica and B. bassiana, five strains for each genetic 
group were examined for the width and the length of conidia 
and phialides (Table S3). The data were analysed with one- 
factor ANOVA for testing difference between species, as well 
as between the genetic groups within each species. 


RESULTS 


Population structure and molecular diversity 


The maximum-likelihood tree inferred from all the 729549 
SNPs showed clear differentiation between B. asiatica and 
B. bassiana as expected, and also intraspecific genetic sub- 
division (Fig. 1a). Two and three clades could be distinguished 
within B. asiatica and B. bassiana, respectively. One of the two 
clades of B. asiatica corresponds to the Clade C as found in 
Khonsanit et al. (2020), including the ex-type strain ARSEF 
4850; the other clade comprises a mix of individuals from the 
Clades A-D from that previous study, therefore named here as 
'Clade Mixed'. The three clades from B. bassiana correspond to 
the Clades A-C established in this species by Khonsanit et al. 
(2020); the Clade C includes the ex-type strain ARSEF 1564. 
The principal component analysis (PCA) also showed clear 
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Fig. 1 Population genomics analyses. a. A circular NJ tree based on whole-genome SNPs, the red asterisks represent the type strains of Beauveria asiatica 
(ARSEF4850) and B. bassiana (ARSEF 1564), thick branches correspond to those supporting 100 96 bootstrap; b. principal component analysis (PCA) for total 
sampling; c. Bayesian clustering for total sampling; d. Bayesian clustering for B. asiatica; e. PCA for B. bassiana; f. PCA for B. asiatica. 
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Table 1 Genetic differentiation (F, below the diagonal), genetic divergence (D,,: above the diagonal) between clades of B. asiatica and B. bassiana, and 


nucleotide diversity (Pi: italic values forming the diagonal) within clades. 


B. asiatica B. bassiana 
Clade C Clade Mixed Clade A Clade B Clade C 
B. asiatica Clade C 0.009 0.0875 0.651 0.675 0.669 
Clade Mixed 0.742 0.036 0.630 0.654 0.648 
B. bassiana Clade A 0.991 0.969 0.003 0.157 0.119 
Clade B 0.970 0.949 0.893 0.031 0.161 
Clade C 0.980 0.958 0.913 0.848 0.018 


difference as expected between B. asiatica and B. bassiana 
as well as between these two species and the others included 
(Fig. 1b). Not enough individuals had been analysed from the 
other species to assess their differentiation from each other. 
The Bayesian clustering within each of the two species recov- 
ered the same genetic clusters as above (Fig. 1c—d) but, for B. 
asiatica, an additional subdivision was observed within Clade 
Mixed (Fig. 1d). The PCAalso confirmed the strong differentia- 
tion between the clades identified in the analyses above (Fig. 
1e-f) as well as the subdivision within B. asiatica Clade Mixed 
(clusters hereafter called mix1 and mix2: Fig. 1f). These sub- 
clusters included mixes of strains from various genetic groups 
from Khonsanit et al. (2020). The various analyses thus consist- 
ently indicated the existence of intraspecific subdivision. The 
only exception was the strain NHJ10436 which clustered with 
B. bassiana Clade B in the Bayesian clustering and the tree, but 
had intermediate coordinates on the PCA (Fig. 1b, e). This is 
probably due to missing data as this strain only carried 14427 
SNPs (1.97 % from total SNPs). This strain was though kept in 
our analyses as it personally interested us for use in biocontrol. 
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The F, between the two species was high (0.875) as expected 
for interspecific differentiation. The nucleotide diversity (Pi) was 
higher in B. bassiana (0.107) than B. asiatica (0.057); this may 
be due to the use of B. bassiana genome as reference to call 
SNPs. The differentiation levels between clades within species 
were also very high, with F... reaching more than 0.7 (Table 1), 
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Fig. 2 A Bayesian phylogenetic tree based on concatenated 1034 genes 
in Beauveria species. The red dots represent the nodes supported by 1.00 
posterior probability. The highlighted clades are supported by Bayesian 
concordance factors with 95 96 confidence interval. 
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Fig. 3 Morphometric analysis. a. Boxplots representing the distribution of 
the length and the width of phialides and conidia between Beauveria asiatica 
and B. bassiana. The red asterisks denote significant difference between 
the two species; b. boxplots representing the distribution of the length and 
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Fig. 4 Statistical analyses of secondary metabolite profiles resulting from 
liquid chromatography-mass spectrometry (LC-MS). a. Euclidean distance- 
based NJ tree; b. principal component analysis (PCA) for total sampling 
(top panel) and for the clades within Beauveria asiatica and B. bassiana 
(bottom panel). 


but the F.. between the sub-clusters mix1 and mix2 within Clade 
Mixed of B. asiatica was relatively weak (0.262). The divergence 
estimated by D, was also relatively high between clades, close 
to or higher than 0.1 (Table 1), but that between the sub-clusters 
mix1 and mix2 was much lower (0.041). The nucleotide diversity 
(Pi) was higher in B. asiatica Clade Mixed (0.024 and 0.036 for 
the sub-clusters mix1 and mix2, respectively) than in Clade C 
(0.009). In B. bassiana, Clade B had the highest Pi (0.031), 
followed by Clade C (0.018) and Clade A (0.003). 


Phylogenetic network and linkage disequilibrium 


The Neighbour-Net network inferred from p-distances between 
the strains showed clear separation between clades within each 
species, with minimal reticulation between them, confirming 
lack of gene flow among clades (Fig. S1). The sub-cluster mix1 
was nested within the sub-cluster mix2 for B. asiatica Clade 
Mixed. Reticulations were observed within B. asiatica Clade C 
and Clade Mixed. For B. bassiana, reticulation events could be 
observed only within Clade B. In all the clades within species, 
the r? dropped to rapidly reach a plateau except B. bassiana 
Clade A which sustained a longer LD with r? dropping under 
half the maximum around 300 Kbp (Fig. S1). These results are 
consistent with the absence of reticulation in the B. bassiana 
Clade A, for which the nucleotide diversity was also the lowest 
(0.002: Table 1). 


Phylogenetics species recognition 


The species tree inferred from 1034 concatenated genes re- 
covered strongly supported monophyletic clades corresponding 
to the genetic subdivisions within B. asiatica (Clade C and 
Clade Mixed) and B. bassiana (Clades A-C) (Fig. 2). The sub- 
cluster mix 1 within B. asiatica Clade Mixed was recovered as 
a well-supported monophyletic clade but not the sub-cluster 
mix 2 which did not form a monophyletic clade, further sup- 
porting the view that mix1 is a lineage having originated from 
mix2. To evaluate whether these genetic subdivisions could 
be considered as distinct species using phylogenetics species 
criteria, we conducted Bayesian concordance analyses between 
genes having no missing data for B. asiatica with its outgroup 
B. brongniartii (ARSEF617) (887 genes), and for B. bassiana 
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Species 


asiatica C (B. asiatica s.s.) 


asiatica Mixed (B. thailandica) | 
bassiana A (B. namnaoensis) 


bassiana B (B. neobassiana) 
bassiana C (B. bassiana s.s.) 


Presence/Absence 


WAbsence Presence 


NHJ13051 


Cluster! tepene 
Others Cluster! nps 


Clustert7 tipksaps 
Cluster21 tepene 
Cluster23 other 
a Cluster2a tipks 
Cluster25 terpene 
Cluster26 lepene 
mm Cluster27 nps 
Cluster28 tipks 
Cluster29 nps 
Clusters tipks 
= Cluster30 terpene 
m Clusters NPS 
Cluster32 terpene 
Cluster33 NIPS 
Cluster34 tipks-nps 
Cluster35 t1pks 
Cluster36 nips 
[] Glustera7 indolenms 
Cluster3g other 
Cluster39 NIPS 
Cluster4 t3pks-tlpks 
Cluster40 nips 
Cluster41 tepenenps 
Clusters tipks 
Cluster lepene 
Cluster2 PS 


Cluster10 tipks-nmps 
Clusteri2 tipks 
Clustert3 tipks-nps 
Clusteri4 nps 
Cluster15 other 
Clusterí other 
Clusteri8 nps 
Clusteri9 tipks 
Cluster20 other 
Cluster22 nps 


Fig. 5 The presence/absence of secondary metabolites gene clusters among Beauveria species. The isolates were ordered from left to right according to a 


binary distance-based clustering. 


with its outgroup B. kipukae (ARSEF7032) (100 genes). We 
found that the Clades A-C within B. bassiana were highly 
concordant throughout the genome with high Bayesian con- 
cordance factors (BCF), NHJ10436 being well anchored in 
Clade B. For B. asiatica, Clade C had a high BCF (» 0.8) but 
Clade Mixed had a low BCF (0.201) (Fig. 2). 


Morphological variation 


The average measurements of the length and the width of 
phialides and conidia are reported per strain in Table S3. The 
average length and width of phialides and conidia of B. asia- 
tica s.lat., B. bassiana s.lat. and all clades inside each of 
them can be found in Table S4. The length of phialides and 
conidia were significantly different between B. asiatica and 
B. bassiana (one-factor ANOVA; phialide length: F = 16.858, 
p-value = 4.324e; conidia length: F = 79.566, p-value = 
6.309e) but the other traits were not significantly different 
between the two species (phialide width: F = 0.199, p-value = 
0.659; conidia width: F = 79.566, p-value = 6.309e°°) (Fig.3a). 
Between Clade C and Clade Mixed of B. asiatica, only the length 
of conidia was significantly different (F = 7.219, p-value = 0.02, 
Fig. 3b). For B. bassiana, none of the examined traits was signi- 
ficantly different between clades (Table S4) (Fig. 3b). 


Chemotaxonomy 


Metabolite profiles obtained with LC-MS are shown in Fig. S2— 
S3. Principal component analyses (PCA) and the Euclidean 
distance-based tree showed a relatively clear separation be- 
tween B. asiatica and B. bassiana with some exceptions (Fig. 4); 
the two strains MY11641 and NHJ13734 of B. asiatica showed 
a chemical similarity to B. bassiana. The other species included 
in the chemical analyses were not all well discriminated from 
B. asiatica and B. bassiana. On the NJ tree, B. caledonica, 
B. pseudobassiana, B. sungii and B. vermiconia clustered within 
B. asiatica (Fig. 4a). Based on the PCA, B. amorpha, B. mala- 
wiensis, B. pseudobassiana, B. varraoe could be well distin- 
guished from both B. asiatica and B. bassiana. In contrast, 
B. gryllotalpidicola and B. kipukae appeared nested within 


B. bassiana while B. caledonica, B. sungii and B. vermiconia 
were found very close to B. asiatica (Fig. 4b, top panel). The 
intraspecific clades were not well separated based on these 
methods. 


The distribution pattern of secondary metabolite production cor- 
responded broadly to the presence/absence pattern of secon- 
dary metabolite gene clusters (SMGCs) among the species 
(Fig. 5). There were 41 SMGCs inferred from the reference 
genome (Table S4, S5), of which 15 were conserved among 
all the isolates. Twenty-five Beauveria isolates possess all the 
SMGCs most of which belong to B. bassiana s.lat. The three 
clades within B. bassiana s.lat. did not show distinct clustering 
patterns with most of the strains maintaining all SMGCs, except 
Cluster 37 (Indole-Nrps) that was randomly lost among some 
strains (Fig. 5). The isolates of B. asiatica s.lat. grouped together 
with a few SMGCs deletions shared between them — Cluster1 
(Terpene), Cluster11 (NRPS), Cluster41 (Terpene-NRPS). 
No distinction could be observed between the two clades of 
B. asiatica. The other species included tended to group together 
outside the two species complexes but some of them were 
scattered among B. asiatica (Fig. 5). 


Taxonomy 


Based on genomic and morphological results above, B. asi- 
atica Clade Mixed is proposed as B. thailandica; B. bassiana 
Clade A and Clade B are proposed as B. namnaoensis and 
B. neobassiana, respectively. 


Beauveria namnaoensis Khons., Kobmoo & Luangsa-ard, 
sp. nov. — MycoBank MB 838940; Fig. 6 


Etymology. Name derived from the location where the type specimen 
was found, Nam Nao National Park, Nam Nao District. 


Holotypus. THAILAND, Phetchabun Province, Nam Nao District, Nam Nao 
National Park, Headquarter Nature Trail, N16.74° E101.57°, on adult of Xylo- 
copa latipes (Hymenoptera, Apidae), 4 July 2012, A. Khonsanit, K. Tasanathai, 
P. Srikitikulchai, S. Mongkolsamrit & W. Noisripoom (holotype BBH 36158; 
ex-holotype strain BCC 64218). 
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Fig. 6 Beauveria namnaoensis. a—d. Fungus on the hosts; e. colony obverse on PDA after 10 d; f. colony reverse on PDA after 10 d; g. colony obverse on 
PDA after 20 d; h. colony reverse on PDA after 20 d; i-j. phialides and conidia; k. conidia. — Scale bars: a, e—h = 10 mm; b-d = 5 mm; i = 10 um; j-k = 5 um. 


Sexual morph — Unknown. 

Asexual morph — Hosts covered with white mycelium, powdery 
when sporulating. Phialides hyaline, solitary, smooth-walled, 
base ampulliform, clavate, globose, mucronate, pyriform, sub- 
globose, subspherical to lageniform (2—)2-3.2(—5) x (1.5—) 
1.9-2.7(—3.5) um (n = 150). Conidia hyaline, smooth-walled, 
globose to subglobose, occasionally obovoid, (1.5—)1.9-2.8 
(—4) x (1.5—)1.7-2.3(-3) um (n = 150). 

Colony characteristics — Colony growing at room tempera- 
ture attaining 29-31 mm diam in 10 d, and 42-43 mm in 20 d. 
Surface mycelium dense, floccose, sporulation starts at 4 d after 
inoculation, white with green-yellow (1D), becoming to green- 
yellow (1C) and powdery while sporulating. Colony reverse 
green-yellow (1A), white at the margin. 

Hosts — Specimens found on blue milkweed beetle (Coleop- 
tera), stink bug (Hemiptera, Pentatomidae), wasp (Hymenop- 
tera), pupa and adult of moth (Lepidoptera) and grasshopper 
(Orthoptera). 


Additional specimens examined. Tua AND, Chanthaburi Province, Khao Soi 
Dao District, Khao Soi Dao Wildlife Sanctuary, N13.10° E102.19°, on stink bug 
(Hemiptera, Pentatomidae), 1 Jan. 2000, B. Thongnuch, J.J. Luangsa-ard, 
K. Tasanathai, P. Srikitikulchai, S. Mongkolsamrit & W. Chaygate (BBH 
14343: BCC 18114); Kanchanaburi Province, Thong Pha Phum District, 
Thung Yai Naresuan Wildlife Sanctuary, Krathon Ruesi Nature Trail, N15.33° 
E98.92°, on moth adult (Lepidoptera), 2 Dec. 2005, B. Thongnuch, K. Tasa- 
nathai, P. Srikitikulchai & W. Chaygate (BBH 15191: BCC 19745); Phet- 
chaburi Province, Kaeng Krachan District, Kaeng Krachan National Park, 
Ban Krang Camp Nature Trail, N12.88? E99.63*, on blue milkweed beetle 
(Coleoptera), 26 Apr. 2006, B. Thongnuch, K. Tasanathai, R. Ridkaew & 
W. Chaygate (BBH 18375: BCC 21293); Nakhon Ratchasima Province, 
Pak Chong District, Khao Yai National Park, Bueng Phai Nature Trail, 
N14.44° E101.37°, on wasp (Hymenoptera), 29 Nov. 2006, B. Thongnuch, 
P. Puyngain, T. Keokene & W. Chaygate (BBH 23082: BCC 23823); Chiang 
Mai Province, Doi Inthanon National Park, Mae Chaem District, Mae Chaem 
Junction (KM.38) Nature Trail, N18.54° E98.52°, on grasshopper (Orthop- 
tera), 27 Nov. 2008, A. Khonsanit, K. Tasanathai & P. Srikitikulchai (BBH 
25297: BCC 34350). 


Notes — Beauveria namnaoensis is hardly distinguishable 
morphologically from other closely related species such as 


Table 2 Morphological comparisons of Beauveria bassiana species complex. 


Species Host/substrate Phialides (um) Conidia (um) Distribution References 

B. bassiana s.str. Beetle adults, cucurbit beetle, 2.5-6 x 3-6 2-3x2-3 Brazil, Canada, Hungary, Rehner et al. (2011), 
mantid, maple leafcutter moth larva, Morocco, Republic of Korea, Khonsanit et al. (2020), 
moth adult, mulberry moth pupa, Thailand, USA This study 
spittle bug, weevils 

B. namnaoensis sp. nov. Blue milkweed beetle, grasshopper, 2-5 x 1.5-3.5 1.5-4 x 1.5-3 Thailand This study 
moth, lepidopteran pupa, stink bug, 
wasp 

B. neobassiana sp. nov. Ant, beetle adults, cicada adults, 2-3.5 x 2-3 2-3 x 1.5-2.5 Thailand This study 


lepidopteran larva, queen ants, 
moth, weevils 
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B. bassiana and B. neobassiana (Table 2 and refer to the results 
on morphological variation). However, this new species is sup- 
ported as distinct based on phylogenetic analyses. 


Beauveria neobassiana Khons., Kobmoo & Luangsa-ard, 
sp. nov. — MycoBank MB 838939; Fig. 7 


Etymology. Name derived from the morphological similarity to Beauveria 
bassiana. 


Holotypus. THAILAND, Nakhon Ratchasima Province, Pak Chong District, 
Khao Yai National Park, trail from Heo Sawat to Heo Sai Waterfall Nature Trail, 
N14.44° E101.37°, on adult cicada (Hemiptera), 20 May 1996, K. Tasanathai, 
N.L. Hywel-Jones, S. Sivichai & S. Thienhirun (holotype BBH 5054; ex- 
holotype strain BCC 1848). 


Sexual morph — Unknown. 

Asexual morph — Hosts covered with white mycelium, powdery 
while sporulating. Phialides hyaline, solitary, smooth-walled, base 
ovoid, subspherical (2—)2.3 -3.1(—3.5) x (2—)2.1-2.9(—3) um 
(n = 150). Conidia hyaline, smooth-walled, globose, occasion- 
ally subglobose, (2-)2-2.6(-3) x (1.5-)1.8-2.2(-2.5) um 
(n = 150). 

Colony characteristics — Colony on PDA at room tempera- 
ture attaining a diam of 20 -21 mm in 10 d, 40-41 mm in 20 d, 
surface mycelium dense, floccose. Sporulation starts 5 d after 
inoculation, white with green-yellow (2D), powdery when sporu- 
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lating. Colony reverse yellow (6C), becoming to yellow-orange 
(20A) and white at the margin. 

Hosts — Specimens found on adult beetles, weevils (Co- 
leoptera), adult cicadas (Hemiptera), ants (Hymenoptera) and 
larvae and adult of moth (Lepidoptera). 


Additional specimens examined. THAILAND, Kanchanaburi Province, Sang- 
khla Buri District, Khao Laem National Park, Sum Nuk Buhd To, N15.02° 
E98.60°, on adult beetle (Coleoptera), 21 June 1995, N.L. Hywel-Jones, 
R. Nasit & S. Sivichai (BBH 4614: BCC 1446); Phetchaburi Province, Kaeng 
Krachan District, Kaeng Krachan National Park, Khlong 1 Nature Trail, 
N12.88° E99.63°, on moth (Lepidoptera), J.J. Luangsa-ard, K. Tasanathai, 
N.L. Hywel-Jones, P. Srikitikulchai, R. Nasit, S. Mongkolsamrit & W. Chay- 
gate (BBH 14301); Chiang Mai Province, Mae Chaem District, Doi Inthanon 
National Park, N18.54° E98.52°; on ant (Hymenoptera), 1 Jan. 2000, 
B. Thongnuch, K. Tasanathai, P. Srikitikulchai, R. Ridkaew, S. Mongkolsamrit 
& W. Chaygate (BBH 14373: BCC 18124, BBH 14374: BCC 18171); Phetch- 
abun Province, Nam Nao District, Nam Nao National Park, trail to nature study 
area, N16.74° E101.57°, on weevil (Coleoptera, Curculionoidea), 23 May 
2000, K. Tasanathai, P. Lutthisungneon, R. Nasit & S. Sivichai (BBH 7666: 
BCC 2660); Tak Province, Umphang District, Umphang Wildlife Sanctuary, Pi 
Tu Kro Waterfall Nature Trail, N15.92° E98.76°, on queen ant (Hymenoptera), 
26 June 2001, A. Khonsanit, B. Thongnuch, J.J. Luangsa-ard, K. Tasanathai, 
P. Srikitikulchai, S. Mongkolsamrit & W. Chaygate (BBH 23919: BCC 31619); 
Chiang Mai Province, Fang District, Doi Phahompok National Park, Doi 
Phahompok Nature Trail, N20.00° E99.14°, on queen ant (Hymenoptera), 
5 Feb. 2006, B. Thongnuch, K. Tasanathai, S. Mongkolsamrit & W. Chaygate 
(BBH 16599: BCC 20197); Nakhon Ratchasima Province, Pak Chong District, 
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Fig.7 Beauveria neobassiana. a—h. Fungus on the hosts; i. colony obverse on PDA after 10 d; j. colony reverse on PDA after 10 d; k. colony obverse on PDA 
after 20 d; |. colony reverse on PDA after 20 d; m-n. phialides and conidia; o. conidia. — Scale bars: a—b, d-h = 5 mm;c,i-l = 10 mm; m = 10 um; n-o = 5 um. 
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Khao Yai National Park, Kong Kaeo Waterfall Nature Trail, N14.44° E101.37°, 
on adult beetle (Coleoptera), 12 June 2007, C. Chuaseeharonnachai, S. Sivi- 
chai & S. Mongkolsamrit (BBH 24543: BCC 25950); Chanthaburi Province, 
Soi Dao District, Khao Soi Dao Wildlife Sanctuary, Headquarter Nature Trail, 
N13.10° E102.19°, on adult cicada (Hemiptera), 22 July 2008, A. Khonsanit, 
J.J. Luangsa-ard, K. Tasanathai, P. Srikitikulchai & S. Mongkolsamrit (BBH 
24402: BCC 30545), on adult beetle (Coleoptera), 22 July 2008, A. Khon- 
sanit, J.J. Luangsa-ard, K. Tasanathai, P. Srikitikulchai & S. Mongkolsamrit 
(BBH 23856: BCC 31604); Chiang Rai Province, Chiang Khong District, 
Phlu Kaeng Waterfall Nature Trail, N20.26° E100.39°, on Lepidoptera larva 
(Lepidoptera), 17 Jan. 2009, A. Khonsanit, K. Tasanathai, P. Srikitikulchai 
& S. Mongkolsamrit (BBH 26635: BCC 35935); Nakhon Ratchasima Pro- 
vince, Pak Chong District, Khao Yai National Park, Mo Sing To Nature Trail, 
N14.44° E101.37°, on adult beetle (Coleoptera), 20 July 2009, K. Tasanathai, 
P. Srikitikulchai, R. Ridkaew, S. Mongkolsamrit & T. Chohmee (BBH 26358: 
BCC37372); ibid., 16 Aug. 2009, K. Tasanathai, P. Srikitikulchai, R. Ridkaew, 
S. Mongkolsamrit & T. Chohmee (BBH 26797: BCC 37938); Chanthaburi 
Province, Pong Nam Ron District, Suchin's Orchard, N12.92° E102.39°, on 
adult beetle (Coleoptera), 16 Dec. 2014, A. Khonsanit, D. Thanakitpipattana, 
N. Wiriyathanawudhiwong & W. Noisripoom (BBH 39718: BCC 76579). 


Notes — Beauveria neobassiana is hardly distinguishable 
morphologically from other closely related species, such as 
B. bassiana and B. namnaoensis (Table 2 and refer to the 
results on morphological variation). However, this new species 
is based on phylogenetic studies. 


Beauveria thailandica Khons., Kobmoo & Luangsa-ard, 
sp. nov. — MycoBank MB 838941; Fig. 8 


Etymology. The epithet refers to the locality where the type specimen 
was found, Thailand. 


Holotype. THAILAND, Chanthaburi Province, Soi Dao District, Khao Soi Dao 
Wildlife Sanctuary, N13.10° E102.19°, on Coleoptera larva (Coleoptera), 
1 Jan. 2000, N.L. Hywel-Jones (holotype BBH 13831; ex-type culture BCC 
16585; TBRC 8350). 


Sexual morph — Stromata arising from the posterior part of Co- 
leoptera larva, cylindrical with rounded apices, yellow (8B-8A) 
when fresh and yellow-orange (19B -19A) when dry, 38-50 x 0.8 - 
2.5 mm. Perithecia semi-immersed, ovoid, (420—)469.6 - 539 
(—580) x (150—)182.2- 258.5 (- 290) um (n = 30). Asci hyaline, 
cylindrical, capitate, (180—)215.4 -299.5(-335) x (3-)3.2- 
4.7(-6) um (n = 50). Asci-caps hyaline, hemispherical, (2—) 
2.2-3.1(-3) x (3.5—-)3.6-4.1(-4) um (n = 30). Ascospores 
hyaline, filiform, multiseptate, 241-320 x 1 um (n = 8), breaking 
into 64 part-spores. Part-spores hyaline, cylindrical with trun- 
cated end, (4—)5.2-12.9(-23) x 1 um (n = 50). 

Asexual morph — Hosts covered with white mycelium, powdery 
when sporulating. Phialides hyaline, solitary, smooth-walled, 
base ampulliform, mucronate, rostrate, subspherical to lageni- 
form (2—)2.5-3.7(—5) x (1.5—)2-3(-4) um (n = 150). Conidia 
hyaline, smooth-walled, obovoid, ovoid, occasionally globose 
to subglobose (2.5—)3.1-4.3(—6) x (1.5—)1.6-2.3(-3) um 
(n = 150). 

Colony characteristics — Colony growth at room tempera- 
ture attaining a diam of 21-22 mm in 10 d, 30-32 mm in 20 d. 
Surface mycelium dense, convex to the agar surface, floccose, 
cottony, sporulation starts at 10 d after inoculation, white to 
green-yellow (1C). Colony reverse green-yellow (1C), white at 
the margin. 

Hosts — Specimens found on ladybugs, longhorn beetle 
(Coleoptera), earwigs (Dermaptera), adult cicada, leafhoppers, 
stink bug (Hemiptera), queen ant and wasps (Hymenoptera). 


Additional specimens examined. THAILAND, Phetchabun Province, Nam 
Nao District, Nam Nao National Park, Lum Nam Cheun Nature Trail, N16.74? 
E101.57°, on stink bug (Hemiptera, Podopidae), 11 Oct. 1994. N.L. Hywel- 
Jones, R. Nasit & S. Sivichai (BBH 4362: BCC 1442); Phetchaburi Province, 
Kaeng Krachan District, Kaeng Krachan National Park, KM. 15 on road to Tor 
Tip Waterfall, N12.88° E99.63°, on adult cicada (Hemiptera), 23 May 1995, 
N.L. Hywel-Jones (BBH 4541: BCC 1654); ibid., on adult beetle (Coleoptera), 
23 May 1995, N.L. Hywel-Jones (BBH 4542: BCC 1655); ibid., on earwigs 


(Dermaptera, Chelisochidae), 25 May 1995, N.L. Hywel-Jones (BBH 4563: 
BCC 2044); Kanchanaburi Province, Sangkhla Buri District, Khao Laem 
National Park, Sum Nuk Buhd To, N15.02° E98.60°, on adult beetle (Coleop- 
tera), 21 June 1995, N.L. Hywel-Jones, R. Nasit & S. Sivichai (BBH4561: 
BCC1665); Chanthaburi Province, Khao Soi Dao District, Knao Soi Dao 
Wildlife Sanctuary, Behind office on the Nature Trail, N13.10° E102.19°, on 
Lepidoptera larva (Lepidoptera), 20 June 1996, K. Tasanathai, R. Nasit & 
S. Sivichai (BBH 5182: BCC 1906); Phetchabun Province, Nam Nao District, 
Nam Nao National Park, behind office to bungalow 1 Nature Trail, N16.74° 
E101.57°, on Coleoptera larva (Coleoptera), 20 Aug. 1996, K. Tasanathai, 
R. Nasit & S. Sivichai (BBH 5227: BCC 2086); Phetchaburi Province, Kaeng 
Krachan District, Kaeng Krachan National Park, KM. 15 on road to Tor Tip 
Waterfall, N12.88° E99.63°, on weevil bug (Coleoptera, Curculionidae), 
2 June 2000, N.L. Hywel-Jones (BBH 7769: BCC 2676); Mae Hong Son 
Province, Mae Hong Son Road marker KM. 5.7, N19.30° E97.97°, on earwigs 
(Dermaptera, Chelisochidae), 30 June 2002, N.L. Hywel-Jones & R. Nasit, 
(BBH 16766: BCC 12907); Surat Thani Province, Phanom District, Khao 
Sok National Park, Sip Et Shin Waterfall Nature Trail, N8.99° E98.63°, on 
adult beetle (Coleoptera), 30 Sept. 2003, K. Tasanathai, N.L. Hywel-Jones & 
S. Sivichai (BBH 5353: BCC 2120); Phetchabun Province, Nam Nao District, 
Nam Nao National Park, Headquarter Nature Trail, N16.74° E101.57°, on 
adult leafhoppers (Hemiptera, Cicadellidae), 30 June 2004, B. Thongnuch, 
K. Tasanathai & W. Chaygate (BBH 10140: BCC 16183); Tak Province, 
Umphang District, Umphang Wildlife Sanctuary, Pi Tu Kro Waterfall (Preto Lo 
Su), N15.92° E98.76°, on adult longhorn beetle (Coleoptera, Cerambycidae), 
26 June 2008, A. Khonsanit, J.J. Luangsa-ard, K. Tasanathai, P. Srikitikul- 
chai & S. Mongkolsamrit (BBH 23918: BCC 31618); Nakhon Ratchasima 
Province, Pak Chong District, Khao Yai National Park, Mo Sing To Nature 
Trail, N14.7116666666667 E101.421666666667, sexual morph emerging 
between head and thorax of adult click beetle (Coleoptera, Elateridae), 
18 June 2009, N.L. Hywel-Jones, K. Tasanathai, R. Ridkaew, S. Mongkol- 
samrit & T. Chohmee (BBH 38825: BCC 36657); Chanthaburi Province, 
Khao Soi Dao District, Khao Soi DaoWildlife Sanctuary, N13.10° E102.19°, 
on longhorn beetle (Coleoptera), 21 Aug. 2009, B. Thongnuch, J.J. Luang- 
sa-ard, K. Tasanathai, P. Srikitikulchai, S. Mongkolsamrit & W. Chay- 
gate (BBH 14344: BCC 18115); Nakhon Ratchasima Province, Pak Chong 
District, Khao Yai National Park, Mo Sing To Nature Trail, N14.44° E101.37°, 
on Coleoptera larva (Coleoptera), 27 Sept. 2011, A. Khonsanit, K. Sansatcha- 
non, K. Tasanathai, S. Mongkolsamrit & W. Noisripoom (BBH 32175: BCC 
49762); Chiang Mai Province, Chiang Dao District, Chiang Dao Wildlife 
Sanctuary, Chiang Dao Wildlife Research Station, N19.39° E98.84°, on 
Lepidoptera larva (Lepidoptera), 5 Oct. 2012, A. Khonsanit, K. Tasanathai, 
P. Srikitikulchai, R. Promharn & W. Noisripoom (BBH 38847: BCC 56283); 
Chiang Mai Province, Chiang Dao District, Ban Hua Thung Community For- 
est, N19.39° E98.84°, on adult wasp (Hymenoptera), 31 Oct. 2014, A. Khon- 
sanit, D. Thanakitpipattana, K. Tasanathai, P. Srikitikulchai, S. Wongkanoun 
& W. Noisripoom (BBH 40622: BCC 76509). 


Notes — Phylogenetically, B. thailandica is closely related 
to B. asiatica. It occurs on a wide range of insect hosts. Mor- 
phologically, B. thailandica occasionally has smaller perithecia, 
shorter asci and ascospores as well as longer part-spores than 
B. asiatica (Table 3). Detailed statistical analyses showed a 
significant difference in terms of conidial length and a notable 
higher variation for the length and the width of phialides and 
conidia. 


DISCUSSION 


Beauveria is a genus with many cryptic species. There has been 
no comprehensive comparison between genetic, morphological 
and chemical data. In this study, our objectives were to assess 
the species status of genetic groups found within B. asiatica 
s.lat. and B. bassiana s.lat. based on population genomics, 
morphological and chemical data. 


Phylogenetics species criterion 


The intraspecific genetic groups as found in these two species 
complexes by Khonsanit et al. (2020) were recovered by popu- 
lation genomics data in this study. The main question is whether 
they should be considered as distinct species. The evolutionary 
species concept (ESC) defines a species as 'a single lineage 
of ancestor-descendant populations which maintain its identity 
from other such lineages, and which has its own evolutionary 
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Fig. 8 Beauveria thailandica. a. Stromata on host (Coleoptera larva); b—g. fungus on the hosts; h. fertile head; i. perithecia; j. asci; k. asci with asci caps; 
|. ascospore; m. part-spores; n. colony obverse on PDA in 10 d; o. colony reverse on PDA in 10 d; p. colony obverse on PDA in 20 d; q. colony reverse on 
PDA in 20 d; r-s. phialides and conidia; t. conidia. — Scale bars: a, n—q = 10 mm; b-g = 5mm; h = 1 mm; i = 200 um; j, | = 50 um; k, r-t = 5 um; m = 10 um. 
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Table 3 Morphological comparisons of Beauveria asiatica species complex. 


Distribution References 


Conidia 
(um) 


Phialides 
(um) 


Part-spores 
(um) 


Ascospores 
(um) 


Asci-caps 
(um) 


Perithecia Asci 


(um) 


Stromata 
(mm) 


Host 


Species 


(um) 


Rehner et al. (2011), 


China, 


2.5-6.3x3-6 2.5-4.3 x 2-3 


228-455x1 3-12x1 


2-3x3-5 


300-630 x 130-330 268-513 x 3-5 


beetle adult, 27-50 x 1.5-5 


cicada adult, 


B. asiatica s.str. 


Khonsanit et al. (2020) 


Republic of Korea, 


Thailand 


Click beetle, dark 


black chafer beetle, 
grasshopper eggs, 


longhorn beetles 


This study 


Thailand 


2-5x1.5-4 2.5-6 x 1.5-3 


241-320x1  4-23x1 


2-3x3.5-4 


420—580 x 150-290 180-335 x 3-6 


38-46 x 0.8-2.5 


B. thailandica sp. nov. cicada adult, 


click beetle, 


earwigs, ladybugs, 


leafhoppers, longhorn 
beetle, queen ant, 
stink bug, wasps 


tendencies and historical fate' (Wiley 1978). This concept is not 
debated but is not an operational species delimitation criterion 
(De Queiroz 2007, Giraud et al. 2008). For fungi, the so-called 
phylogenetic species concept (PSC) is an operational concept 
that has been predominant and proposed as the most consist- 
ent to ESC (Taylor et al. 2000). According to the theoretical 
concept related to the most common ancestor and descending 
progenies, PSC is naturally associated to the monophyly and 
consistency of clades across genes, which became the basis to 
phylogenetic species recognition (PSR) (Cracraft 1983, Mishler 
& Brandon 1987, Mishler & Theriot 2000). According to this 
criterion, B. asiatica Clade C, Clade Mixed and B. bassiana 
Clade A-C should be considered as distinct species. However, 
monophyly-based species recognition received also some criti- 
cisms; monophyly-based species recognition cannot apply to 
asexual species because, without recombination, all lineages 
are monophyletic and with consistent gene genealogies. Beau- 
veria had actually long been thought to be an asexual genus 
(Basalmo-Crivelli 1835). Our analyses of linkage disequilibrium 
(LD) decay suggested substantial recombination within some 
of these clades, which is consistent to observations of sexual 
morphs for several Beauveria species, including B. bassiana 
(Rehner et al. 2011, Khonsanit et al. 2020). The phylogenetic 
networks which take into account reticulated evolution still 
showed clear separation between them. 


However, PSC can be arbitrary in regard of the limit of where 
to put species boundaries as any kind of molecular markers 
with polymorphisms would allow grouping individuals into 
monophyletic clades based on their respective alleles. PSC 
was thus proposed to rely on the concordance between distinct 
gene genealogies (Avis & Ball 1990, Baum & Shaw 1995), re- 
sulting in the genealogical concordance phylogenetic species 
recognition (GCPSR) criterion for which the limit of species 
is placed where there is a transition from the concordance 
among branches connecting different species to the conflict 
between branches within species due to intraspecific recom- 
bination (Taylor et al. 2000). Some fungal species were in fact 
erroneously defined based on only concatenated analysis of 
multi-gene phylogenies without concordance between mark- 
ers and any corroboration from morphological, chemical and 
ecological data (Liu et al. 2016). Our results showed that the 
different genetic clades actually had high concordance factors 
except for B. asiatica Clade Mixed. Except for this latest case, 
the clades as revealed by whole-genome data are thus sup- 
ported as distinct phylogenetic species. 


Insights from the chemotaxonomy and morphology 


Chemotaxonomy is an approach of classifying and identifying 
microorganisms based upon the similarities and differences 
in biochemical compositions. In filamentous fungi, chemo- 
taxonomy in the broadest sense usually involves secondary 
metabolites which have been extensively studied particularly 
in Ascomycota. Penicillium was the first asexual genus to be 
chemotaxonomically examined (Frisvad 1981). Therein, thin- 
layer chromatographic profiles of mycotoxins in combination 
with classical taxonomy allowed the recognition of four new 
groups within Penicillium subg. Penicillium. Later, several 
genera including Alternaria (Andersen et al. 2008), Aspergillus 
(Frisvad & Samson 1990, Samson et al. 2004, Frisvad & Larsen 
2015), Talaromyces (Frisvad et al. 1990), Fusarium (Thrane & 
Hansen 1995, Schmidt et al. 2004, Zain 2010), Stachybotrys 
(Andersen et al. 2003), Trichoderma (Thrane et al. 2001) and 
many more, were subjected to chemotaxonomic examinations 
showing highly species-specific metabolic profiles. The chemo- 
taxonomic concept has also been proven to be successfully ap- 
plied in polyphasic taxonomy within xylarialean fungi which led 
to the discovery of potential chemotaxonomic markers (Stadler 
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etal. 2001a, b, 2003, 2014, Stadler & Hellwig 2005, Kuhnert et 
al. 2017, Kuephadungphan et al. 2021). For instance, sporothric 
acid, isosporothric acid and dihydroisosporothric acid appeared 
to be specific to Hypoxylon monticulosum (Surup et al. 2014), 
viridistratin A-C to Annulohypoxylon viridistratum (Becker et al. 
2020), Minutellins A-D to Annulohypoxylon minutellum (Kuh- 
nert et al. 2017) and lenormandin A-G to Hypoxylon jaklitschii 
and Hypoxylon lenormandii (Kuhnert et al. 2015). 


Considering Beauveria, chemotaxonomy was adopted and ap- 
plied more than three decades ago (Mugnai et al. 1989, Bridge 
et al. 1990), revealing chemical profiles corresponding more 
or less to species based on morphological criteria, but were 
heterogeneous and overlapping between populations within 
species (e.g., B. bassiana). Recently, Berestetskiy et al. (2018) 
demonstrated that the chromatographic profiles could be used 
to distinguish B. bassiana from B. pseudobassiana while Yin 
et al. (2020) expressed doubts on using chemotaxonomy for 
species identification in Beauveria. The latter study showed 
that, although the two species appeared to produce structurally 
different congeners of the cyclodepsipeptide beauveriolides, 
they also produced similar compounds found in fungi from other 
genera. This is not surprising since several fungal species, not 
only Beauveria, can have abundant compounds in common. 
In our study, we found a discriminant pattern of secondary 
metabolites distribution between B. asiatica and B. bassiana, 
but not between the clades within each species. This pattern 
reflected a conservation of all secondary metabolite gene 
clusters shared between most of the B. bassiana s.lat. strains 
while some were lost in B. asiatica s.lat. and other species. 
However, between the clades within each species complex, 
no distinctive discriminatory pattern was found. 


Itis still too early to conclude that the chemotaxonomic concept 
has failed to discriminate these clades apart as the secondary 
metabolites were only obtained under a single growth condition. 
The YMG medium used in our study has been proven to allow 
expression of many species-specific metabolites in various 
ascomycetes and is even used to reveal the existence of cryp- 
tic species in certain genera, i.e., Xylaria (Stadler et al. 2003, 
Kuephadungphan et al. 2021). The metabolic profiles derived 
from YMG seemed to allow discrimination only between some 
Beauveria species in our study. As culture media composition 
has influences on secondary metabolite production, other media 
should be evaluated. For now, it can thus be concluded that 
the metabolic profiles generated under these conditions could 
only be used to distinguish B. asiatica from other Beauveria 
species with some exceptions. 


Previous studies have shown that the size and the shape of 
conidia can be used, to some extent, to discriminate between 
species but these characters were largely overlapping across 
species of Beauveria (Rehner & Buckley 2005, Imoulan et 
al. 2017, Khonsanit et al. 2020). Our study revealed that the 
morphological differences among clades within B. asiatica and 
B. bassiana were much less evident than the genomic data. 
Nevertheless, the morphology, particularly the conidial length, 
was significantly different between B. asiatica Clade C and 
Clade Mixed (conidial length). Altogether, our findings support 
the idea that the two Beauveria species are composed of sibling 
species with similar phenotypes, i.e., cryptic species. 


Taxonomy 


Given the evidence discussed above, we propose B. bassiana 
Clade A and Clade B as new species, named as B. namnaoen- 
sis and B. neobassiana, while Clade C which includes the ex- 
type strain ARSEF 1564 can be considered as B. bassiana s.str. 
Regarding B. asiatica s.lat., the Clade Mixed is proposed here 
as a new species, B. thailandica, and the Clade C containing 
the ex-type strain ARSEF4850 remains B. asiatica s.str. The 


description and photographic materials for the novel species are 
given above in the results. These species are largely overlap- 
ping in their distribution, e.g., sympatric species (Fig. S4), but 
are genetically well isolated with limited recombination between 
them. Despite lack of distinctive distribution pattern of second- 
ary metabolites, the genomic, morphological and ecological 
data supported that they are distinct species. 


Future studies 


Investigating the sexual compatibility between closely related 
species within each of the two species complexes should be 
conducted in future to assess whether they could fit the criterion 
of reproductive isolation for biological species (De Queiroz 
2005). Reproduction mode and mating systems are extremely 
variable in fungi, from pure clonality to outcrossing sexual 
reproduction, with significant implications in genetic diversity 
and adaptability (Billiard et al. 2011, Taylor et al. 2015). Some 
studies in fungal pathogens have already shown that closely 
related cryptic species could be intersterile (e.g., Calonectria 
spp: Lombard et al. 2010, Li et al. 2020) or sexually compatible 
with post-mating isolation (e.g., Microbotryum spp.: De Vienne 
et al. 2009). Beauveria was initially proposed to be an asexual 
genus, i.e., being present naturally only under asexual form 
(Basalmo-Crivelli 1835). Sexual reproduction structures could, 
however, occasionally be observed in natural habitats (Sanjuan 
et al. 2014, Khonsanit et al. 2020). Our study furthermore sug- 
gested some extent of genetic recombination and sexual re- 
production in all the newly proposed species. 


CONCLUSIONS 


Our work showed that whole-genome sequence data could 
provide strong support for intra-species genetic groups and be 
powerful for delimiting species. During the last decade, geno- 
mics data have greatly contributed to the elucidation of species 
complexes in many organisms (Chan et al. 2017, Cerca et al. 
2021) including fungi (Sepülveda et al. 2017, Kobmoo et al. 
2019, Matute & Sepülveda 2019). For hypocrealean entomo- 
pathogenic fungi, such approaches should be particularly use- 
ful as there are likely many cryptic species in different genera 
(Blackwellomyces and Cordyceps: Mongkolsamrit et al. 2020b, 
Isaria: Mongkolsamrit et al. 2018, Metarhizium: Mongkolsam- 
rit et al. 2020a, Ophiocordyceps: Kobmoo et al. 2012, 2019, 
Araüjo et al. 2018, Khonsanit et al. 2019). Many taxonomic 
studies were based on a few samples with a limited number 
of markers which has likely led to the underestimation of the 
diversity in this fungal group. Combining whole-genome data, 
thorough morphological examination with statistical analyses 
and chemical profiling will contribute to establish a solid basis 
for species discovery in fungi. 
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Supplementary material 


Fig. S1 Analyses of recombination footprint for: a. Beauveria asiatica; b. B. 
bassiana. The left panel represents phylogenetic networks based on SNPs 
from respective species. The right panel represents the analyses of linkage 
disequilibrium (LD) decay. 

Fig. S2 Liquid chromatography-mass spectrometry (LC-MS) profiles of 
mycelia extracts from Beauveria asiatica and B. bassiana species complexes. 
Fig. S3 Liquid chromatography-mass spectrometry (LC-MS) profiles of 
mycelia extracts from Beauveria spp. included for comparison with B. asiatica 
and B. bassiana. 


Fig. S4 Distribution map of Beauveria asiatica s.lat. (left panel) and B. bassi- 
ana s.lat. (right panel). 


Table S1 List of Beauveria strains used in the study. The species were at- 
tributed according to previous studies or routine identification in our laboratory. 
Specimens of B. asiatica were proposed to be a novel species, B. thailandica, 
while some of B. bassiana were proposed as new species, B. neobassiana 
or B. namnaoensis sp. nov., based on phylogenomics results. 


Table S2 Fermentation periods for the production of secondary metabolites 
of Beauveria strains used to study chemotaxonomy. 


Table S3 Average width and length (um) of conidia and phialides per strain 
for the strains used in statistical analyses of morphological traits. 

Table S4 The results of statistical analyses on morphological traits. a. 
The average length and width of phialides and conidia among Beauveria 
species; b. ANOVA tests of difference between clades within B. asiatica 
and B. bassiana. 


Table S5 Characteristics of secondary metabolite gene clusters inferred on 
the genome reference of Beauveria bassiana ARSEF8028. 


